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Abstract 

We study the metastable states in Ising spin models with orthogonal interac- 
tion matrices. We focus on three realizations of this model, the random case and 
two non-random cases, i.e. the fully-frustrated model on an infinite dimensional 
hypercube and the so-called sine-model. We use the mean-field (or tap) equations 
which we derive by resuming the high-temperature expansion of the Gibbs free 
energy. In some special non-random cases, we can find the absolute minimum of 
the free energy. For the random case we compute the average number of solutions 
to the tap equations. We find that the configurational entropy (or complexity) 
is extensive in the range T RSB < T < T M . Finally we present an apparently unre- 
lated replica calculation which reproduces the analytical expression for the total 
number of TAP solutions. 
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The aim of this paper is to study the mean field equations (the tap equations) for the 
local magnetization for the fully frustrated Ising model on an hypercube, or equivalently 
on a single cell of an hypercubic lattice in the limit of infinite dimensions. 

The model is very interesting and it has been widely studied. It belongs to the wide 
class of models which have a non-random Hamiltonian but they behave in a similar 
(or identical) way to other random models. These models can be studied using the 
usual techniques for random systems (e.g. the replica method) and in this way one can 
obtain the correct results (with maybe the exception of the equilibrium behavior at low 
temperature) [I], 

In this model (as in the other models of the same class) there are at least two 
transitions: 

• The dynamical transition at which the correlation time diverges. At this temper- 
ature static (i.e. equal time quantities) are fully regular. 

• An equilibrium transition at which the replica symmetry is broken. Below this 
temperature there are many equilibrium states available to the system 

If we cool an infinite system below the dynamical transition, its energy does not 
go to the equilibrium energy and the system remains trapped in a metastable state of 
higher energy. 

Generally speaking one can associate to each stable or metastable state a solution 
of the tap equation. On the contrary the inverse relation does not hold. There are 
many solutions of the tap equations which do not correspond to stable of metastable 
states. According to the folklore a tap solution correspond to a metastable state only 
if it is separated by other solutions by a barrier whose height diverges with the volume; 
it seems most of the solutions are separated by other solutions by finite barriers. 

The relation between metastable states and the exponentially large persistency time 
in a given metastable state from one side, and the properties of the solutions of the TAP 
equations has never fully understood. The aim of this paper contribute to this direction 
by computing some of the properties of the solutions of the tap equation in a model 
which has metastable states. 

An other interesting property of the model (absent in its random correspondent) 
is the presence of a very low energy state, which cannot be reached with continuity 
coming from the high temperature phase. At a given temperature the system has a first 
transition to a state with very low energy (the crystal state) . The behavior of the system 
in this low temperature phase may be understood with good precision by considering 
the corresponding tap equation. 

There is also an other point on which we would like to call the attention of the reader. 
For the Sherrington-Kirkpatrick (sk) model it has already noticed that the computation 
of the total number of solutions of the tap equations turns out to be equivalent to a 
apparently unrelated replica computation. More precisely it was found that: 



AT ~ lim Z. 



m ■ 




2 



where prefactors have been neglected and only the exponentially large terms have been 
taken into account Z m is the partition function computed by breaking the replica sym- 
metry into two groups of replicas, one with m elements, the other with n — m terms. 

The relation is quite surprising, because the r.h.s. may be evaluated without having 
to write down the very TAP equation. 

Apparently there is no known explanation for this phenomenon. Our contribution 
to this point is to observe that the previous relation holds also in this case, where the 
form of the tap equation is much more complex. It is quite like (as also suggested 
by Cugliandolo and Kurchan), that there should be an isomorphism of two algebraic 
structure which explain this equality, but we have not done progresses in this direction. 

It may be interesting to recall that a similar phenomenon happens in ferromagnetic 
systems with a random temperature or magnetic field. Let us consider the case of 
random magnetic field. Here we are interested find the probability of having more than 
one solution to the stochastic differential equation 

- A0(x) +m 2 (f)(x) +g<f)(x) 3 = h(x), (2) 

which plays the role of the tap equation in this model. 

In the replica approach || £|, || it was found that the existence of many solutions is 
related to the presence of non trivial saddle points, where the field </) a (x) for a = l,n 
(n = 0) is of the form 

4> a {x) = f(x) for 1 < a < m 

<fi a (x) = fix) H for m < a < n, (3) 

m 

and one consider the limit where m — > oo. 

Also in this case the replica symmetry is broken by dividing the replicas into two 
groups, one with m elements, the other with n — m terms, and taking the limit m — > oo. 
This phenomenon seems to be quite widespread. 

The plan of the paper is the following. 

In section two we recall the definition of three different models, which have the same 
high temperature expansion: the random orthogonal model (ROM) where the coupling 
matrix is a random orthogonal matrix, the sine model, where the elements of the coupling 
matrix can be written as the sine of an appropriate expression, and the fully frustrated 
models on a single hypercubic cell (ffm). These three models are identical in the high 
temperature expansion, because they have very similar coupling matrices. 

In section three we use the high temperature expansion to derive the TAP equations. 

In section four we study the property of the solution of the tap equation correspond- 
ing the lowest lying state (of zero energy) which exists in the sine model and likely exists 
in the fully frustrated hypercube (the so called crystal state). 

In section five we compute the average number of tap solution for ROM. We study 
their properties as function of their free energy both at T = and T ^ 0. We also 
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find a relation among the properties of the tap equations and the marginality condition 
for the dynamical transition and show that the replica symmetric free energy can be 
written as the sum of contributions from a large number of metastable states above the 
RSB transition. 

In section six we present the computation of the total number of tap solutions with 
the two replica group method we have described before. 

Finally in the appendix we present some technical details needed in the computations 
of section five. 

1 The Model 

In what follows we will consider the model defined by the Hamiltonian 

H = -§53 J ijWj> ( 4 ) 

ij 

where {<7j} is a set of N Ising spin variables (crj = ±1) and Jij is an N x N symmetric 
orthogonal matrix with large connectivity (j>1). To lighten the notation, the matrix 
Jij will be taken to have zeros on the diagonal^. 

1.1 Random Orthogonal Model 

In the random orthogonal model (ROM) the coupling matrix is chosen at random in 
the set of orthogonal symmetric matrices. The probability distribution (or integration 
measure) is defined by writing J = ODO -1 with D a diagonal matrix composed of 
±1 and O a generic orthogonal matrix (not necessarily symmetric) whose probability 
distribution is defined by the Haar measure on the orthogonal group. We will make use 
of the identity f| 




which holds in the large N limit when A is a symmetric matrix of finite rank and where 
G(x) is given by 

1 ,/i _|_ 4 T 2 i i i 1 

G(x) = -~ log( Vi + ^ +i ) + - A VT+4^ - \. (6) 

Notice that integration over matrices chosen from a Gaussian distribution (sk model) 
also yields Eq. flD but with G{x) = x 2 /4. 

1 Strictly speaking a generic matrix J in ROM and the one for the sine model have non-zero diagonal 
elements, nevertheless the matrix obtained by setting those elements to zero will also be orthogonal in 
the large TV limit. 
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1.2 Sine Model 



The coupling matrix for the sine model is given by 

Ja = = sm . (7) 

The matrix J is obviously symmetric and its orthogonality follows from orthonormality 
relations among simple harmonics. This model was introduced in as a simple Hamil- 
tonian that admits a complex ground state for special values of N. More precisely, if 
2N + 1 is prime and iV is odd then the Legendre sequence 

Oj=f (mod2iV + l), (8) 

which consists of ±1 is a ground state configuration of the sine model with energy 
density —1/2. This fact is by no means obvious; the interested reader is referred to 
for details. 

Monte Carlo simulations have shown that the thermodynamical properties of the 
sine model are the same as that of ROM. The only difference is that for the sine model 
with 2N + 1 prime there exists a low-lying state — never seen for large N when cooling 
from high temperature — in which the system remains (at low temperature) if it is put 
there by hand. One explanation of the similarity between the sine model and ROM 
is that Eq. (|7|) can be viewed as a (bad) pseudo-random generator, and therefore the 
couplings are, for all practical purpose, random. 



1.3 Fully- Frustrated Model 

Frustration in an Ising-like system (with Jjj = ±1) is defined by the product of the 
couplings over a given plaquette 0. If this product is —1 the plaquette is said to 
be frustrated. One can construct a coupling matrix on a d dimensional simple cubic 
lattice such that every plaquette is frustrated. This construction is not unique but all 
realizations are gauge equivalent so the thermodynamical properties of the model are 
unique. For a single hypercubic cell, the fully-frustrated condition imposes || 

JikJkj = d5 i:j . (9) 

k 

Therefore if we divide the couplings by yd the matrix J will be orthogonal and sym- 
metric. The thermodynamical limit will be taken by letting d go to infinity, which will 
also insure that the coordination number is large. 

Note that unlike SK on a cubic lattice, for the fully-frustrated model (ffm) the 
d —>■ oo limit of a single hypercube is not completely equivalent to that of the lattice. 
The reason is that for SK the distribution of the eigenvalues of the coupling matrix (the 
celebrated Wigner semi-circle law) is the same both on the single cell and on the lattice. 
In the fully-frustrated case, the matrix J has only two eigenvalues { — 1, 1} for the unit 
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hypercube, while it has a continuous distribution for the lattice, in the range [—2, 2] 
with rms value v^2 ||. Heuristically, in the high-T phase the system is sensitive to all 
eigenvalues of the matrix: as we shall see, both models (lattice and unit cell) have the 
same free energy up to a rescaling of (3 by a factor of \/2. On the other hand, in the 
low-T phase, it is the largest eigenvalues of J that matter: for example the ground state 
of the lattice has twice the energy density of that of the unit cell. It is this difference 
between 2 and \[2 that makes it impossible to go from a single cell to the lattice with 
a simple rescaling of the temperature: their qualitative behavior should nevertheless be 
the same. 

The ground state configuration of FFM is not known for d > 8. It is conjectured 
that, with the \fd normalization, the ground state energy density goes to —1/2 for large 
d and that it is equal to —1/2 when d is a perfect square. As with the sine model, 
Monte Carlo simulations of FFM are well described by ROM; furthermore they show 



aging effects, which strengthens our belief that FFM is a glassy system [10 



2 Mean-Field Free Energy 

2.1 Standard High- Temperature Expansion 

Our goal is to derive the tap equations for the orthogonal model using the high- 
temperature expansion of the magnetization-dependent (Gibbs) free energy. Before 
doing so, let us try to understand the standard (m; = 0) high-temperature expansion. 
We can write 

e"^ = 2Met 1 / 2 (/3J) /n^expj-^^^V^ + E^gcosh^ + ^j, (10) 

where we have introduced a site-dependent magnetic field hi for later convenience; in 
this section we set it to zero. We can view this partition function as a theory of N fields 
in zero dimension. In this language the high-temperature free energy F(/3) is the sum 
of connected diagrams whose vertex factors are the coefficients of the Taylor series for 
log cosh a;, with the propagator /3Ji~. The first few terms are 




where Ju = excludes 'petals' — terms involving J2 Ja- The only terms that survive 
the large N limit are the so-called 'cactus' diagrams — trees made out of loops joined 
at the vertices (see Fig. 1). The sum of Jy's for cactus diagrams is equal to N if the 
lengths of all loops are even (even-cacti) and zero otherwise. To show this, we take a 
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cactus diagram and contract each loop to a point starting from the outermost loops, the 
contraction is done using the orthogonality relation 

'^JikJkj = (12) 

k 

In other words, the high-temperature free energy is given by N times the sum over 
the combinatorial factors and powers of f3 of all even-cactus diagrams. In particular this 
free energy is independent of the particular choice of the orthogonal matrix J. 

There is more than one way to calculate this sum. Since it is independent of the 
matrix J it must be equal to the (annealed) average over such matrices, recovering the 
result of Ref. || 

-PF((3) = N\og2 + NG(P), (13) 

with G(x) given by Eq. @. 

More directly, we notice that the desired diagrams are precisely those summed by 
the high-temperature series of the spherical model. In the spherical model the Ising 
variables Oi are replaced by continuous variables Si constrained to have J2i Sf = N. 
Fixing this last constraint with a Lagrange multiplier, the sum over all configurations 
is then a Gaussian integral which can be done exactly knowing the eigenvalues of the 
matrix J, i.e.,1 and —1 in equal proportions. After eliminating the Lagrange multiplier, 
we find 

- f3F sph ((3) = j (log27r + 1) + NG(J3), (14) 

with the same G(x). The first terms in Eq. (|T3| ) and Eq. (|H| ) are the volume of phase 
space of the two respective models. A similar analysis was done in Ref. |J for ffm on 
the full lattice. They also found Eq. ( pID but with G(j3) replaced with G(y20). 

For comparison, the high-temperature series for SK is easily done; the only surviving 
diagram is the first one in Eq. (O) leaving us with 

-/3F SK (/3) = iVlog2 + iV^, (15) 
where we have used that Jfj is self-averaging. 



2.2 Gibbs Free Energy 

The high-T free energy we have just calculated does not tell us anything about the low- 
temperature phase. It does not have any singularity for positive (3 which might signal a 
phase transition. What we need to compute is the magnetization-dependent free energy 
(often called the Gibbs free energy). In SK this quantity is exactly the TAP free energy 
which is still valid in the low-temperature phase. The program is simple, we need to 
expand Eq. ( |T0D in powers of j3 and hi and perform a Legendre transform, passing from 
the variables hi to using 

dF 

§{(3,mi) = F(P,hi) -^Krrii with mi = —. (16) 
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Calculating the first few terms in (3 by hand, we find 



m<) = -!£{(l + TOi )log[f(l + mi 



+ (1 -rm) log -fin)]) 



+^ £ Jamfinj + jE4( 1 " - m )) + °(^)- ( 17 ) 



The first term is the entropy of Ising spins constraint to have magnetization {mj, 
the second one is minus the 'naive' mean-field energy, and the third corresponds to the 



Onsager reaction term. In Ref. |TT| the tap equations for SK were re-derived in this way, 
the higher order terms being negligible in that case. The analysis of the previous section 
(rrii = 0) tells us that for a model with an orthogonal matrix, the high-temperature 
series will include an infinite number of non-negligible terms. We must therefore find a 
systematic way of calculating all the terms in the Gibbs free energy. 

In standard field theory, the Legendre transform can be formulated in diagrammat- 
ical terms. The Gibbs free energy, also called the effective potential in particle physics, 
is given by the sum over connected one-particle-irreducible diagrams. This is not true 
in our case. The reason is that the external field hi does not appear as a linear source 
term in Eq. (TTH), instead it appears inside the potential. A shift in <pi would not simplify 
matters since it would introduce a term quadratic in hi sufficient to render invalid the 
usual 1PI derivation. The present authors don't know of any diagrammatic expansion 
of the Gibbs free energy for Ising-like system. 

It is however possible to perform such an expansion algebraically and to give it a 
diagrammatic interpretation [Q. The weak point of this method is that the vertex 



weight and the combinatorial factors cannot be calculated systematically. Nevertheless, 
one can see that the cactus diagrams^], will appear in the Gibbs free energy with the 
same weight but with an extra factor of (1 — mj.) for each vertex i% 
have 



13 . Therefore we 



~\ Y, {(! + m *) l0 § [K 1 + m *)l l0 § 



'1 - rrii 



'J 



ij 

P 4 



rri: 



/3 4 



^(l-m^)4(l-m|)j| fc (l-m 2 fc ) 

ijk 



+ x £(! - rr$)Ml ~ m^Ml - ml)J kl {\ - mf)J u + 0(/3 6 ). (18) 

8 ijkl 

We now argue that in the large iV limit, all those terms (except for the 'entropic' and 
'energetic' one) are self-averaging, that is, the error introduced by substituting ml by 
q = iV -1 Y!,i m i vanishes in the thermodynamical limit. Equivalently in SK the Onsager 

2 What we call cactus diagrams are called 'loop diagrams' in Ref. [Q. The difference in terminology 
arises from the fact that these authors consider restricted sums (sums over different indices with com- 
binatorial factors) while we consider unrestricted sums. The vertex weights coming from the expansion 
of log cosh a; are precisely those necessary to go from one type of sum to the other. 
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term in the free energy is often replaced by /3 2 (1 — g) 2 /4. We can now resum Eq. (IS), 
by noticing that we are now summing over the same diagrams than in the m = case. 
We recover the reaction term of Eq. (|13| ) with (3 replaced by /3(1 — g), 



Z ij 



|(1 + m 



+ (1 - 772;) log |(1 - mi) } 



(19) 



where the function G(x) is again given by Eq. (|j). 

The same result can be obtained in a more heuristic way. We saw in the previous 
section that the spherical model reproduces the high-temperature expansion of the model 
with Ising spins except for the entropic term. It is quite plausible that the two models 
will have the same Gibbs free energy, except once again for the entropic term. This 
argument has been used (the other way around) to deduced the right tap equations for 
the spherical p-spin model JTJJ. The general Gibbs free energy for a spherical model 
with quadratic interaction is given by 



sph 



logdet 1/2 (A - (3 J) - ^(log27r + A) + W m *( X k 



3 i 



(20) 



with A determined by the saddle point equation. Specializing to the case where the 
eigenvalues of J are {+1,-1} in equal proportions and after eliminating A we find 



/3$ 



sph 



-y(log[27r(l 



?)] + x ) - § S J ij m t m j 

ij 



NG(P(l-q)), 



(21) 



with q = N~ Y^i m i- The first term is the entropy of spherical spins constrained to have 
magnetization mj. As expected, it is the only term that differs from Eq. ([T9]). Note that 
in this case we did not need to substitute 1 — m| with its average value (1 — q). 

Armed with the free energy (|19"1), we can now write down the 'mean-field' or TAP 
equations for our model. They are given by the partial derivatives of $ with respect to 
the m^'s: 



tanh 



rru 



2/3G'(J3(l-q))m i -/3 y £j i 



yTUj 



0. 



(22) 



We note finally that if we substitute in Eqs. (19, |22D G(x) = x 2 /A we recover the 
standard result of []T3| for SK just as the same substitution in our integration formula (^j) 
recovers the integration over a Gaussian distribution. In Sect. f|, this fact will provide 



us with an easy way to compare our formulae against those of Ref. [16 



3 Minimum of the Free Energy 

The orthogonality of J imposes a simple bound on the energy of any spin configuration. 
Indeed, a configuration vector {crj} can be decomposed into its projection {s^} on the 
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eigenspace with eigenvalue +1 and {s { } on the complement with eigenvalue 
we have 



E(a) 



1 







2 






( 


S+ 




S 


1 



-1. Then 



(23) 



In other words, the energy is bounded below by —N/2 and this bound is reached iff a 
is an eigenvector of J with eigenvalue +1. We shall argue later that for a generic large 
orthogonal matrix the existence of an eigenvector whose every component is either +1 
or —1 is highly improbable. Nevertheless, we can construct such a matrix, for example 
the sine model for odd N with 2N + 1 a prime number. Whether or not ffm admits 
such a ground state for special values of d and/or in the limit d — > oo is still an open 
question. Unfortunately, neither an e = —1/2 spin configuration for d > 9 has been 
found nor has a proof that it cannot exist. 

For the remaining part of this section, we will consider a model that admits an 
e = — 1/2 ground state. The tap equations (^) admit in this case a solution of the 
form 

™* = Vq^i, (24) 

where the {e{\ are ±1 and form an eigenvector of J with eigenvalue +1. With this 
ansatz Eq. (|22]) reduces to 



q = tanh 2 < j3y/q 



1 + 



1 + 4/? 2 (l - f?) 1 



(25) 



where we have used Eq. 
given by 



2/3(1 - q) 

to compute G'(x). This solution has specific free energy 



(3f 



1 + 



log |(i + 



+ 



log |(1 " y/q) 



?))■ (26) 



For a given value of q, this is the lowest free energy solution. In fact, the energetic 
term reaches its lower bound from orthogonality, the entropy is maximum (at fixed q) 
when all the local magnetizations are equal in magnitude and the reaction term only 
depends on q. The absolute minimum of the free energy must therefore be of this form. 
Eq. (P5|) always admits a solution with q = 0. At low temperature it will also admit 
a solution with non-zero q. Numerically one finds that below T = 0.4000 a solution 
with q = 0.92 appears but with higher free energy that the paramagnetic solution. At 
T = 0.178 (q = 0.99995) this solution becomes the true minimum of the free energy. 

that the system undergoes a first order transition at 



One might conclude as in [13 
= 0.178. 



T = 0.178. This transition is not seen in Monte Carlo simulations of ffm or the sine 
model when cooling down from high temperature. Besides Eq. (p4"|), there are many 
other solutions to the mean-field equations. As we shall we in the next section, for a 



3 This paper contains many numbers obtained from numerical analysis, they are noted with equal 
signs and every digit is significant. 
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generic interaction matrix and below a certain critical temperature, their number grows 
exponentially with the size of the system. It is those solutions and the large free energy 
barriers between them that prevent the system from finding the true minimum of the 
free energy and gives it instead a glassy behavior. 

When the ground state configuration is known — as in the sine model for special N— 
it is however possible to see this 'crystalline state' in Monte Carlo simulations. One has 
to start in the ground state at T = and continuously increase the temperature. Fig. 
10 of Ref. || shows the result of such a simulation (sine model N = 44 and N = 806). 
Calculation of the internal energy for our solution reproduces exactly those curves. 
Note that for the larger value of iV the system stays trapped beyond the point where 
the paramagnetic solution becomes stable (T = 0.712 on their scale) this is a clear sign 
that barriers between different metastable states are very high. 

The solution can be transposed directly to the ffm on the lattice by taking care 
of the appropriate factors of 2 and y2 recovering the result of Ref. |L3| . One important 



point, though, is that this solution is only valid if the ground state energy density is 
— 1/2 (—1 for the full lattice), otherwise {e^} is not an eigenvector of J and Eq. ( |22"D 
does not reduce to Eq. (p5|). 



4 Number of TAP Solutions 



4.1 General Result 

In this section we will compute the average number of metastable states for ROM. We will 



follow as much as possible the steps and the notation of Ref. [|16| where the corresponding 
result for SK was first derive. Recall our tap equations: 



Ti = ttmh- 1 m i + 2[3G'(p(l-q))m i - J < 



i-jTUj 







(27) 



m?1 



with corresponding free energy written as a sum of single site terms: 
/ = ((3N)~ l [~ log 2 - G(/3(l - q)) - /3qG'(j3(l - q)) + \m t tanh" 1 m t + \ log(l 

(28) 

We write the number of solutions as the integral over all possible values of m« fixing 
with ^-functions that the {rrii} form a solution of Eq. ([27]) with free energy /. Using a 
Fourier representation of the 5-function, we obtain 



N 



9 i 1 , f io ° d\ f io ° du f ioc -n ( dxA f 1 . , , r , T/ , 



-A m l + u H /( m ») + H Xi9( m i) - P XI J ij ( x i m i + X i m i) 

i i i i<j 



|det A |(29) 
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where 



.4 



ij 



drrij 



1 _ m 2)-i + 2 /3G"(/3(l - q)j\ % - /3J, 



'J- 



(30) 



dropping a term in rriimj/N (see appendix). 

We can now proceed to average Eq. (|29|) over the random couplings J^-. We should 
really be averaging logjV s , the extensive quantity, not Af s . To do so we would need to 
introduced replicas. The replica symmetric computation would be tedious but straight- 
forward leaving us not with 7 unknown parameters (as will be the case below) but on 
the order of 14. On the other hand, the direct average will be sufficient to provide us a 
clear picture of the metastable states. 

There are two terms that depend on the { Jij}. It is shown in appendix that the two 
can be averaged independently. Using Eq. (^) to average the term involving Xi and rrii, 
we find 



exp 



-j3 E Jij {xirrij + Xjmi 

i<j 



exp <^ NTiG 



x (g) m + m g) x 
N 



exp {iV [G(v + + G(v - Vw)] } , (31) 



where 



N t-r 1 



and w 



N 4- 



(32) 



The second equality in Eq. ( |3"TD follows from an analysis of the eigenvalues of the matrix 
x <S> m + m (g) x. This matrix has only two non-zero eigenvalues corresponding to the 
two terms: v + y/w and v — \/w. 

To compute the average the determinant we will need to drop the absolute value. 
This corresponds to weighting each solution with the sign of the determinant of its 
Hessian matrix. Formally, we would be computing a topological invariant (from Morse 
theory) which as little to do with our original goal. Nevertheless, the calculation with- 
out the absolute value gives sensible results (here and in Refs. KB], |T7[]) and connects 



smoothly with the zero temperature results where the calculation can be done without 
this pathology. This problem is discussed in more detail in Refs. JTB|, [0| . 

We introduce a set of anti-commuting (Grassman) variables {9i,9i} to express the 
determinant as an exponential, 



(det A) = ^ / riOdO exp 

' dOdO exp \ <kdA + NTrG 



y"l cbiOiOi — fizL Jij{0i9j + 6j6i 

* (ij) 

0®Q-Q®d 



0- 



N 



(33) 
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If we expand TrG[. . .] in a Taylor series, we can use the following propriety of Grassman 
variables 



hOi 2 + 9i 2 9 il ){9i 2 9 iz + 9 i3 9i 2 ) . . . ( y 9 in 9 h + 9 h 9 in ) — —2 9 lk 9 ik , (34) 

k=l 



which allows us to resum the series and find 

(det A) = I d0dO exp j ^ - 2NG 



—— \ 9 Si 



(35) 



Finally we fix r = (P/N)J2i9i9i using a Lagrange multiplier R, which allows us to 
perform the integration over the Grassman variables, leaving us with 



/co rioo A T) 

dr — J] (a, + PR) exp{N[-rR - 2G{r)}}. 
-oo J— ioo 27ri 



(36) 



We now collect results form Eqs. (p9|), (|3T|) and (|36|), introduce Lagrange multipliers V 
and to impose Eqs. (P^D, perform the Xj integration, drop multiplicative prefactors 
and finally set V" = 2G'(/3(1 - q)) - A, i2 = B/{3 - 2G'(f3(l - g)) and r = (3b, to obtain 

(jV s (/)) = saddle exp {iV [-Ag -uf-bB-vA- wW/q + + v 7 ^) 

+ G?(u - v 7 ^) - 2C(/36) + 2(/36 - v)G'((3(l - g)) + log/] } , (37) 



where 
/ = 



1 dm 
-l 



2/3vW VI 



i? exp 



tanh 1 m — (3 Am 



+ Am 2 + uf(m) 



(38) 



We have indicated by the term saddle that the right-hand side of Eq. ([37]) is to be 
extremized with respect to the nine variables: A, q, u, b, B, v, A, w, W. By partial differ- 
entiation with respect to those variables, we obtain the saddle point equations. They 
admit a solution with B = which implies 6 = 1 — g; following Ref. (HJ we adopt this 
solution. This choice will also lead us to the correct T = result. 

As advertised earlier, the result for SK can be obtained from Eq. ([J7|) by setting 
G(x) = x 2 /4. In this case, three variables can be eliminated using the saddle point 



equations to give back the same expression as in Ref. [16 . 



4.2 T = Result 

Before looking at the numerical solution for the saddle point of Eq. ([37D , let us study 
its T = limit. In this limit, q — > 1, A decouples, and the five other parameters have 
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a finite value. The integral in Eq. (38) can then be done analytically using a change of 



variable m = 1 — exp(— /3y) and keeping only terms with a finite limit when (3 — > oo, 

lim / = exp(u 2 W/A - Am/2) erfc ( uW ~ ) . (39) 



Shifting A and doing the integral over w, we obtain 

(Af a (f)) = saddle exp [n [-vA - wW + (f + vW + A/2) 2 /W + G(v + y/w) 

+G(v - <Jw) + log erfc(- A/2VW)] } . (40) 

The same result can be obtained directly by counting the average number of spin con- 
figurations {a,} satisfying 

Vi CTj ^2 JijVj > aIld f = Jij a i a 3- ( 41 ) 

j^i i<j 

This computation can be done in rather straightforward way using Heaviside step func- 
tions instead of ^-functions. This approach is free of the problem of the absolute value 
of the determinant. 

The saddle point equations obtained from Eq. ( f4"0"D were solved numerically to give 
the shape of the distribution of one-flip stable configurations as a function of energy 
density (Fig. 2 full-line). Note that since Eq. ( |4*0"D involves the function G(x) for poten- 
tially complex arguments (the complex conjugate term assures us that the final result 
will always be real) one has to choose a sign convention for the square-root terms in the 
definition of G(x) (c.f. Eq. (||)). It was chosen such that the real part of the square-root 
is always of the same sign as v in G(v + ^/w); this way the saddle point never has to 
cross any cut in the complex plane. 

As one follows the saddle point to the rightmost part of the curve (e > —0.36), the 
argument of the error-function in Eq. ([!(]) goes through infinity after which all four 
parameters become complex. It is not clear if the contour of the integral computed by 
the saddle point method can be deformed to have this complex saddle point has its main 
contribution. This issue was not investigated any further for it is of marginal interest. 

In addition, an approximate enumeration of the one-flip stable configurations for 
a relatively small (N = 48) ROM was done. Random configurations were generated 
and then cooled at zero temperature, the resulting one-flip stable configurations were 
compared with previously stored one and were stored if not previously obtained. After 
4 x 10 6 iterations, the low-energy configurations were each found about a hundred times 
while only a few high-energy ones remained unfound. The resulting distribution was 
binned and its logarithm (divided by 48) is shown as the data points on Fig. 2. An 
overall normalization constant (not predicted by the saddle point method) was added 
to the simulation data to make them fall on the predicted curve. Points on the x- 
axis correspond to energy bins where no metastable states were found. The agreement 
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between the theory and this simulation is very good. Similar data (not shown) for 
the sine-model also gave very good agreement with the ROM prediction. For ffm the 
Diophantine constraints (integer spin and integer field) are such that for small d only a 
few values of the energy are allowed for a metastable state (e.g. 2 for d = 5 and 3 for 
d = 6). Those constraints might not be so important for large d, but for d > 6 (N > 64) 
an exhaustive search becomes impossible. Whether the distribution of metastable states 
in FFM follows that of ROM remains an open question. 

The point where the curve shown in Fig. 3 intersects the x-axis on the left, cor- 
responds to the minimum energy (-Ensoi) at which there are on average exponentially 
many configurations. It is almost equal to the approximate ground-state energy (-Ea- S t ep ) 
from the one-step replica breaking solution of this model ||. Careful analysis of the 
equations leading to these two quantity done using arbitrary precision arithmetics shows 
that while they are not equal they differ by less than 30 parts in a billion. Precisely, 
-Ensoi = -0.484119415 and -Ei_ stC p = -0.484119428. The true average ground state 
energy density for ROM is probably well approximated by either of these two quantities. 



By taking partial derivatives of Eq. fl37|) , we can write the number of metastable states 
in terms of the solution of six coupled non-linear equations, three of which containing 
definite integrals. A numerical solution of these equations was done by following the 
T = solution of the previous section to finite T. At fixed T, varying u traces out the 
N^ 1 log(AQ vs / curve. Typical curves (for T = 0.134 and T = 0.065) are shown on 
Fig. 2 (dash line and dotted line); results for other temperature are quite similar. The 
abrupt stop on the right-hand side of these curves is due to the divergence of the definite 
integral in Eq. (|37|) , it is possible that the saddle point can be analytically continued to 
complex values giving a smooth curve all the way to zero as in the T = case. This 
issue has not been investigated since the interesting points lie on the low free energy 
end of the curve. 

Setting u = in Eq. fl3~7|) gives the total number of tap solutions. It gives a positive 
value for iV -1 log(AQ starting from 0.2854 at T = and decreasing all the way to 
at T = 0.32. Above this temperature, the q = saddle point with log(AQ = is the 
correct one. The full results are plotted in Fig. 3 (full line). 

In the thermodynamics of the model nothing special happens at the temperature 
T = 0.32 at which exponentially many tap solutions appear. These metastable states 
have too high free energy to contribute to the partition function. Instead, if we write 



4.3 Finite T Results 



« max e 



/ 



NA(f) 



(42) 



where a labels tap solutions and with 



A(f) = N-Hog(Af s (f))-Pf, 



(43) 
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we find that the solutions which contributes to the free energy are those for which A(f) 
is maximum. They are given by setting u = —/3 in Eq. (|37|) . With this substitution, the 
saddle point equations admit a simple solution of the form: 



2„2 



P 2 q 



A = G'{(3) - G'{(3{1 - q)) 



W — — and A 

P 



And the number of contributing tap solutions reduces to 



where 



./ 



<AQ = exp {N (0 f + G(f3) + log 2)} 



0-1 [- log 2 - G((3(l - q)) - (3qG'((3(l - q)) + (3A 
r°° dz 



(45) 



-0A 



2tt 



=e 2 



cosh(J2/3Az) logcosh(A/2/3A2 



(46) 



with A(q) as above and q determined variationally The function / corresponds the free 
energy of the contributing solutions. We notice that the logarithm of the number of 
contributing solutions, often called complexity, can be written as 



d 

log(AQ = -^f3F 1 _ step (/3,m,q) 



(47) 



m=l 



where -Fi_ ste p is the free energy computed with one-step replica symmetry breaking with 
RSB parameter m, q = and extremized over qi = q (Eq. (44) in Ref. 0). Equality 
([|7]) is in fact more general than what is presented in here, for a more direct derivation 
and its physical interpretation see Ref. |2(J. 

The condition for Eq. (fHf) to have a saddle point at a non-zero value of q is equivalent 
to the marginality condition (T < T M = 0.134). At this temperature, one finds q = 0.962 
with iV -1 log(A/" cs ) = 0.158. Monte Carlo simulations have indicated that all three 
orthogonal models considered here (ffm, ROM and the sine model) undergo a dynamical 
glassy transition at a temperature equal or very near T M 0, [UJ. 

The complexity goes to zero at the replica symmetry breaking transition T RSB = 
0.065. For temperatures below T RSB , Eq. ( |45| ) gives an unphysical negative complexity: 
the above saddle point is no longer valid and the contributing solutions are those with 
the smallest free energy (the leftmost points of the curves of Fig. 2). The full curve of 
complexity vs. temperature is plotted on Fig. 3 (dash line). 

We finally notice that we naturally have 



maxA(/)=log2 + G(/3) = -/3/ H 



(4f 



In other words even in the region between T RSB and T M , where a large number of solutions 
contribute to the free energy, the RS free energy is still valid. The decomposition of the 
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Temp. 


Comments 


U.4U 


Metastable crystalline phase (if present) 


U.oz 


nxponeirciai numuei oi iaf solutions 


0.18 


Stability of crystalline solution (if present) 


0.134 


TAP solutions start contributing to the free energy, 'marginality condition' 




and dynamical phase transition 


0.065 


Complexity goes to zero and replica symmetry breaking 



Table 1: Interesting temperatures for the orthogonal model, from text and adapted from 
Ref. to present normalizations. 



RS free energy into contributions from many metastable states below a temperature "T g " 
higher than T RSB has been discussed in the context of many different models PTl 

Above T M , one can still find a saddle point of Eq. fl3~T| ) with u = —ft, corresponding to 
the maximum of A(f), but it is no longer of the simple form ([H]). Numerically one finds 
that for those temperature, maxA(f) < —ftf RS . The metastable states don't contribute 
to the free energy. 

It may seem surprising that the RS solution is still valid while a very large number 
of different states contribute to the free energy. The explanation is that their number 
is so large that if ones chooses randomly two such states the probability of finding the 
same one is zero, therefore the order parameter P(q) only measures the overlap between 
different states which is always equal to go = 0. Only when the complexity ceases to be 
extensive can the P(q) become non-trivial, leading to replica symmetry breaking. 

We now have a simple intuitive scenario for the phase diagram of models such as 
ROM where the order parameter jumps discontinuously at the transition (see Table 1). 
At a relatively high temperature (T = 0.32 in rom) the number of metastable states 
increases dramatically, but their free energy is so high that they influence neither the 
static nor the dynamics. At T M , a large number of these states start contributing to the 
free energy. At this temperature, all static thermodynamical quantities are perfectly 
regular, but the time scales involved with the dynamics diverge and the system is no 
longer able to thermalize. Phenomena such as aging start to appear. At T RSB the entropy 
of the contributing states goes to zero, P{q) becomes non-trivial and replica symmetry 
is broken. 

5 Two-Group model 

In this final section we will divert a bit our analysis of the thermodynamics of the 
orthogonal model to look at a puzzling analytical 'coincidence.' Indeed, by computing 
the partition function for a certain replica-symmetry breaking scheme, we will recover 
the expression for the average number of tap solution. The 'two-group model' was first 
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introduced in Ref. |25j as an attempt to break replica symmetry in sk. It was later 
noticed |26| that the partition function (lim n ^ (^2G) j) * n this framework is not equal 
to unity but instead one has 

(Zh) = I df(K(f))e- n ^ } (49) 



and in particular lim n ^ (-^2G) J gives the average number of tap solutions ((AQj). We 
present here a general derivation of the partition function for the 'two-group model.' We 
will assume that the coupling matrix is chosen from a distribution that obeys a relation 
such as Eq. @ but we will not make use of the explicit form of the function G(x). We 
will recover the result for SK as the special case G(x) = x 2 /4. 

The average partition function for ROM with n replica is given by 

(Z n )j = JdQdA exp {N [TtG((3Q) - §Tr(AQ) + logZ (A)]} , (50) 

with 

Z (A)=^exp|^A ab( xV b l. (51) 

<J a {a<b J 

The two-group ansatz correspond to writing the matrix A and Q as made up of two 
diagonal blocks of size m x m and (n — m) x (n — m) and equal elements outside those 
blocks. We specialize right away to the case n = which will make some of our formulae 
look a bit strange but will make the notation more compact. The following derivation 
can be extended without much difficulties to the n^O case. We define 



and Q = yi H , (52) 

with A aa = and Q aa = 1 and with the elements parameterized as 

Ui = u + t/m + r/m 2 , u 2 = u — t/m + r/m 2 (53) 

qi = q + a/m + 7/m 2 , q 2 = q — a/m + 7/m 2 , (54) 
where m will be taken to infinity. Having defined the two-group ansatz, we can now 



compute the different parts of Eq. (|50|) ; 

Tr(AQ) = 2{at + qr + yu) - 2(au + qt). (55) 

To compute TrG(/3Q) we find the different eigenvalues of the matrix Q and their mul- 
tiplicity, leading us to 

TrG(/3Q) = lim^ (m-l)G ((3(1 - qi ))-(m + l)G ((3(1 -q 2 )) + G((3q + ) + G((3q-) 
= -2G(f3(l-q))-2a(3G'(f3(l-q)) + G(f3q + ) + G(f3q~), (56) 
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where 

q ± = q-a-l± ^2q(^-a) 
The Z (A) term requires more work, 

1 l uS 2 + (± + -1W + 



(57) 



—As 2 



(58) 
(59) 



where 

S = ^u ni S + =^(j a and S- = ^ a a- 

a<m a>m 

If we make three Hubbard- Stratonovich transformations to linearize the quadratic terms 
and do the sum over spins, we find 



Z (A) 



lim e 

m— >oo 



00 dh 



00 mdh+dh„ 



+m 



3 \flnu J - --- ■' - x 
log cosh {h + /i H 



2vrt 





" h 2 






u 




ft 2- 


)_ 2T" 


H 2f 



+ 



r{h\ + ^) 



(60) 



where we have put in evidence the terms with a factor of m. We can now use the saddle 
point method to compute the m —>■ oo limit of the last expression. To leading order in 
m, the saddle point equations for h + and h_ are the same. They do differ by a term of 
order 1/m but this perturbation only contributes to the saddle point at order 1/m (it 
does not contribute at order one because the derivative of the order m term vanishes). 
At the saddle point, one has h + = h_ = h Q which solves 



E(h, h Q ) = tanh(/i + h Q ) 



T 



0. 



(61) 



The order m contributions to the saddle point of h + and /i_ exactly cancel each other 
leaving us with a finite result as m —>■ oo. Including Gaussian fluctuations, we find 



Z a {A) 



dh 



dE 



-oo tv27ra \dh 
Using the implicit function theorem, 



exp 



2rh 2 ; 



dh 



0h o 



dK vi, 



d ho (i - f 



(62) 



(63) 



we can change the integration variable from h to h Q which we rescalefj to find 

2 



Zo(A) 



1 dm 



-i V2 



7VU 



1 



rrr 



exp 



tanh 1 m — tm 
2u 



+ rm 



(64) 



4 Note that we call m the rescaled integration variable h a /t to make the connection with Eq. (f37|), it 
is not to be confused with the replica-symmetry breaking parameter m which as been taken to infinity. 
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Finally if we collect the results from Eqs. (|55| , 56, 34), make the following change of 
variables 



A = t/P , W = u/2p 2 , A = r , w = p 2 {2q(>y - a)} and v = (3{a + l- q), (65) 

and extremize over the six parameters, we recover (jV s ) given by Eq. Q57| ) with u = 
where the saddle point with B = and 6 = 1 — g has been chosen. Had we not set n = 0, 
we would have recovered the r.h.s. of Eq. (|37l) with it = — n/3 and without the — w/ term. 
In other words (Z% G ) is the Laplace transform of (Af s (f)) with —n[3N conjugate to /. 
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Appendix 

We will now proceed to show that the determinant of the matrix A and the term 
~ JijXinij (c.f. Eq. ( p9|) ) can be averaged separately. Had we kept them together, 
introduced Grassman variables and integrate over the couplings, we would have obtained 
a term like this 



exp <^ iVTrG 



'x®m + m<g)x + 0<g)0-0®0 
* N 



(66) 



Therefore we need to show that the cross-terms in the Taylor expansion of TrG[. . .] are 
irrelevant in the large iV limit. A typical cross-term is 

•^11^12 ^2^*3 ^13^14^14 "^16 " " " ^in^il (^^) 



ATn-l ^ 
ly {ik} 

We can consider those terms as a perturbation about the fermionic Gaussian integral 
( |35f) . We need to compute all connected diagrams generated by these new 6 and 9 
vertices, contracting them with the Gaussian propagator (cij + f3R)~ 1 5ij. The point is 
that since the propagator is diagonal, all such contractions will be at most of order 1 
while the logarithm of the unperturbed result (Eq. ([36])) is of order N. The same is true 
of the term we dropped in writing Eq. ((30). 
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0.3 




0.35 



Figure 3: Normalized logarithm of the average total number of TAP solutions for ROM 
as a function of temperature (iV -1 log(ftf s ) vs T) (full line); same curve but for solutions 
that contribute to the free energy (dash line). 
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Figure 2: Logarithm of the density of metastable states divided by N as a function 
of free energy (iV -1 log(A/" s (/)) vs /) at T = (full line); T = 0.134 (dash line) and 
T = 0.065 (dotted line); Simulation at T = for N = 48 ROM (data points). 
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(a) (b) (c) (d) 

Figure 1: (a) a diagram with a petal, (b) and (c) two non-cactus diagrams, (d) a generic 
even-cactus diagram 
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